Fitting Secondary Growth Models

Statistical Examination of Challenge Test Data and Modelling

2024-04-30

Objectives

  • You will learn how to:
    • Import a .csv file to R software
    • Check a data set
    • Select an adequate secondary model
    • Fit a secondary model
    • Examine the fitting results
    • Obtain the confidence and prediction intervals
  • Materials needed:
    • This presentation
    • The R script: Unit1-ChallengeTestsDataAnalysisSGM.R
    • The accompanying video

Data

Data

  • Growth rate values of salmonella in broth
  • The listeria.csv file has three columns (Temp, GR sqrtGR) and 29 observations
    • Temp is the temperature used to measure the growth rate
    • GR is the growth rate
    • sqrtGR is the square root of the growth rate
  • Lets read the listeria.csv with the read.csv() function and save it as data frame (df)
df <- read.csv("../data/listeria.csv", sep=";", header=TRUE)
head(df)
  Temp      GR sqrtGR
1  0.7 0.00000 0.0000
2  3.4 0.00034 0.0185
3  7.6 0.00127 0.0356
4  9.6 0.00166 0.0407
5 11.2 0.00230 0.0480
6 12.8 0.00256 0.0506

Data

  • Lets, check the dataset using the structure (str()) function
str(df)
'data.frame':   28 obs. of  3 variables:
 $ Temp  : num  0.7 3.4 7.6 9.6 11.2 12.8 14.1 15.4 16.9 18.2 ...
 $ GR    : num  0 0.00034 0.00127 0.00166 0.0023 0.00256 0.00311 0.00423 0.00452 0.00587 ...
 $ sqrtGR: num  0 0.0185 0.0356 0.0407 0.048 0.0506 0.0558 0.065 0.0672 0.0766 ...
  • We confirm that the data set is composed of 3 variables (columns)
    • Temp: temperature
    • GR: growth rate
    • sqrtGR: square root of the growth rate

Check the data shape

  • Check the data by plotting the sqrtGR against Temperature
    • We can observe a non-linear relationship
plot(GR ~ Temp, data=df,
     xlab="Temperature (ºC)", ylab="Growth rate")

plot(sqrtGR ~ Temp, data=df,
     xlab="Temperature (ºC)", ylab="Growth rate")

Secondary models

Secondary models

  • Describe the influence of environmental factors on microbial growth rate
  • Can be generally classified as:
    • Arrhenius-type models
    • Belehradek-type models
    • Square-root models
    • Cardinal parameter models
  • When adjusting a secondary model, using directly the untransformed growth rate response variable produces heterocedastic models
  • To obtain an adequate fitting, the growth rate values are often transformed by taking the natural logarithm or the square root

Model fitting tools

Software

  • We use the R software
    • Its free
    • Developed for statistical analysis and plots
    • Functions to fit non-linear models to experimental data
      • nls() function from the```stats package
      • gsl_nls() from gslnls package (GNU Scientific Library)
      • Others
  • Download and install the R software: https://cran.r-project.org/
  • Install the predmicror package: https://fsqanalytics.github.io/predmicror/

Software

  • Install the predmicror package from GitHub repository
devtools::install_github("fsqanalytics/predmicror")
  • Load the predmicror package
library(predmicror)

Secondary growth model

Selected model: Cardinal model for temperature

library(predmicror)
CMTI
function (x, Tmax, Tmin, MUopt, Topt) 
{
    if (!requireNamespace("gslnls", quietly = TRUE)) {
        stop("Package \"gslnls\" must be installed to use this function.", 
            call. = FALSE)
    }
    CMT <- ifelse(x <= Tmin | x >= Tmax, 0, MUopt * (((x - Tmax) * 
        (x - Tmin)^2)/((Topt - Tmin) * ((Topt - Tmin) * (x - 
        Topt) - (Topt - Tmax) * (Topt + Tmin - 2 * x)))))
    result <- sqrt(CMT)
    return(result)
}
<bytecode: 0x57af77b9bab0>
<environment: namespace:predmicror>

Fitting procedure

Starting values

  • To fit non-linear models we need to we need to supply starting values for the model parameters
  • So, lets start by defining the starting values
start.values = list(Tmax=47, Tmin=0.1, MUopt=0.020, Topt=37)

Fit using gsl_nls() function

  • Now we can fit the CMTI model to the experimental data
library(predmicror)
library(gslnls)
fit <- gsl_nls(sqrtGR ~ CMTI(Temp,Tmax,Tmin,MUopt,Topt),
           data=df,
           start =  start.values
           )
fit
Nonlinear regression model
  model: sqrtGR ~ CMTI(Temp, Tmax, Tmin, MUopt, Topt)
   data: df
    Tmax     Tmin    MUopt     Topt 
44.27284  1.38405  0.02139 37.33238 
 residual sum-of-squares: 0.000484

Algorithm: multifit/levenberg-marquardt, (scaling: more, solver: qr)

Number of iterations to convergence: 6 
Achieved convergence tolerance: 3.619e-09

Check fitting results

  • Now we can summarise and inspect the fitting results by using the summary() function
summary(fit)

Formula: sqrtGR ~ CMTI(Temp, Tmax, Tmin, MUopt, Topt)

Parameters:
       Estimate Std. Error t value Pr(>|t|)    
Tmax  4.427e+01  6.121e-02 723.351   <2e-16 ***
Tmin  1.384e+00  5.564e-01   2.488   0.0202 *  
MUopt 2.139e-02  4.531e-04  47.211   <2e-16 ***
Topt  3.733e+01  3.157e-01 118.245   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.004491 on 24 degrees of freedom

Number of iterations to convergence: 6 
Achieved convergence tolerance: 3.619e-09
  • We can extract the model coefficients using the coef() function
coefs <- coef(fit)
coefs
       Tmax        Tmin       MUopt        Topt 
44.27284300  1.38405426  0.02139205 37.33238243 

Parameters confidence intervals

  • Extract the confidence intervals of the parameters using the confint() function
confint(fit)
        2.5 %      97.5 %
1 44.14652164 44.39916436
2  0.23575472  2.53235380
3  0.02045686  0.02232723
4 36.68076914 37.98399571

Plot the fitted values

  • Define a vector with auxiliary temperature data
new.temp=seq(0.7,44, by=0.2)
  • Use the predict() function to compute the prediction interval
fits <- predict(fit,
                newdata = data.frame(Temp=new.temp),
                interval = "prediction", level = 0.95)
str(fits)
 num [1:217, 1:3] 0 0 0 0 0.000574 ...
 - attr(*, "dimnames")=List of 2
  ..$ : NULL
  ..$ : chr [1:3] "fit" "lwr" "upr"
  • Check the fits object
head(fits[, ])
              fit          lwr         upr
[1,] 0.0000000000 -0.009268608 0.009268608
[2,] 0.0000000000 -0.009268608 0.009268608
[3,] 0.0000000000 -0.009268608 0.009268608
[4,] 0.0000000000 -0.009268608 0.009268608
[5,] 0.0005735376 -0.008828908 0.009975983
[6,] 0.0015625183 -0.007836776 0.010961813

Plot confidence interval

  • Create a plot of the original data with the fitted values superimposed
  • Plot the observed data using the plot() function
  • Use the lines() function to add the confidence interval
plot(sqrtGR ~ Temp, data=df, ylim=c(0,.2))
lines(new.temp, fits[, 1], col="blue")
lines(new.temp, fits[, 2], col="red")
lines(new.temp, fits[, 3], col="red")

To start

Shiny app CardinalFit: https://ubarron.shinyapps.io/CardinalFit/

References

Dolan, Kirk D., and Dharmendra K. Mishra. 2013. “Parameter Estimation in Food Science.” Annual Review of Food Science and Technology 4 (1): 401–22. https://doi.org/10.1146/annurev-food-022811-101247.
Gonzales-Barron, Ursula, and Vasco Cadavez. 2019. Handbook of Predictive Microbiology Growth Models Using R. Edited by Ursula Gonzales-Barron and Vasco Cadavez. Bragança, Portugal: Bringráfica Indústrias Gráficas, LDA.
McMeekin, T. A., J. Olley, T. Ross, and D. A. Ratkowsky. 1993. Predictive Microbiology: Theory and Application. Taunton, UK: Research Studies Press.
Ratkowsky, D A, R K Lowry, T A McMeekin, A N Stokes, and R E Chandler. 1983. “Model for Bacterial Culture Growth Rate Throughout the Entire Biokinetic Temperature Range.” Journal of Bacteriology 154 (3): 1222–26. https://jb.asm.org/content/154/3/1222.
Stephen A., Morse, and Meitzner Timothy A. 2012. “Medical Microbiology.” In, edited by G. Brooks, K. C. Carroll, J. Butel, and S. Morse, 26th Edition, 55–65. Jawetz, Melnick, & Adelberg’s. Mcgraw-hill. https://books.google.pt/books?id=UUSXV8B9i9sC.